Streszczenie

Cel badania

Analiza obejmuje eksplorację zbioru danych zawierającego 925 obserwacji dotyczących elektrod grafenowych stosowanych w superkondensatorach. Dane w zbiorze opisują budowę tych elektrod, uzyskane wydajności oraz inne parametry dotyczące testowania działania elektrod. Głównym celem było zidentyfikowanie kluczowych czynników wpływających na pojemność elektrody (Capacitance) oraz zbudowanie modelu predykcyjnego umożliwiającego prognozowanie tej wartości na podstawie użytych materiałow ich parametrów fizykochemicznych i warunków testowych.

Najważniejsze wyniki

Po oczyszczeniu danych uzyskano finalny zbiór zawierający 908 obserwacji i 12 zmiennych, z zerowym udziałem wartości NaN. Dane charakteryzują się wysoką heterogenicznością oraz brakiem normalności rozkładów. Korelacje między zmiennymi są przeważnie słabe. Warunki testowe odgrywają większą rolę w określeniu pojemności niż same właściwości materiałowe elektrod. Obcięcie wartości ekstremalnych dla 5 zmiennych z >10% outlierów poprawiło stabilność modeli.

Model Random Forest osiągnął współczynnik determinacji R² = 0.7, wykazując najlepszą zdolność predykcyjną spośród testowanych algorytmów (SVR: R² = 0.31, Elastic Net: R² = 0.17). Uzyskał on najniższy błąd RMSE = 237.88 F/g oraz MAE = 137.39 F/g, co potwierdza jego przewagę nad metodami liniowymi w obsłudze nieliniowych zależności. Analiza SHAP wykazała, że warunki eksperymentalne (przedział potencjału, typ elektrolitu) mają większy wpływ na pojemność niż cechy strukturalne materiałów (SSA, Ratio ID/IG).

Konfiguracja środowiska i wykorzystane biblioteki

Wykorzystane pakiety R

packages <- c(
  "tidyverse", "dplyr", "readr", "knitr", "kableExtra", "corrplot", "plotly",
  "randomForest", "caret", "psych", "gridExtra",
  "DataExplorer", "fastshap", "DescTools", "glmnet"
)

new_packages <- packages[!(packages %in% installed.packages()[,"Package"])]
if (length(new_packages)) {
  install.packages(new_packages, dependencies = TRUE, quiet = TRUE)
}

load_package <- function(pkg) {
  if (!require(pkg, character.only = TRUE, quietly = TRUE)) {
    stop(paste("Nie udało się załadować pakietu:", pkg))
  }
  TRUE
}
invisible(sapply(packages, load_package))

Tabela 1. Wykorzystane pakiety i ich funkcje

Pakiet Funkcja w projekcie
tidyverse Zintegrowany zestaw narzędzi do manipulacji i wizualizacji danych
dplyr Operacje na ramkach danych: filtrowanie, grupowanie, agregacja
readr Wczytywanie plików CSV z automatycznym wykrywaniem typów
knitr Generowanie dynamicznych raportów HTML
kableExtra Formatowanie i stylizacja tabel
corrplot Wizualizacja macierzy korelacji
plotly Tworzenie interaktywnych wykresów
randomForest Implementacja algorytmu Random Forest
caret Framework do trenowania i walidacji modeli ML
psych Zaawansowane narzędzia statystyczne (describe)
gridExtra Układanie wielu wykresów w siatce
DataExplorer Automatyczna eksploracja danych (plot_missing)
fastshap Obliczanie wartości SHAP dla interpretacji modeli
DescTools Dodatkowe narzędzia statystyczne (CramerV)
glmnet Modele liniowe z regularyzacją (Elastic Net)

Zapewnienie powtarzalności wyników

set.seed(123)

Ustalenie ziarna generatora liczb pseudolosowych zapewnia pełną powtarzalność wyników przy każdym uruchomieniu analizy.

Wczytanie i inspekcja danych

data <- read_csv("data.csv", show_col_types = FALSE)

cat("Wymiary zbioru danych:", nrow(data), "obserwacji i", ncol(data), "zmiennych\n\n")
## Wymiary zbioru danych: 925 obserwacji i 21 zmiennych
cat("Nazwy kolumn:\n")
## Nazwy kolumn:
cat(paste(names(data), collapse = "\n"))
## Ref.
## Limits of Potential Window (V)
## Lower Limit of Potential Window (V)
## Upper Limit of Potential Window (V)
## Potential Window (V)
## Current Density (A/g)
## Capacitance (F/g)
## Specific Surface Area (m^2/g)
## Charge Transfer Resistance (Rct) (ohm)
## Equivalent Series Resistance (Rs) (ohm)
## Electrode Configuration
## Pore Size (nm)
## Pore Volume (cm^3/g)
## Ratio of ID/IG
## N at%
## C at%
## O at%
## Electrolyte Chemical Formula
## Electrolyte Ionic Conductivity
## Electrolyte Concentration (M)
## Cell Configuration (three/two electrode system)
data %>% 
  head(10) %>% 
  kable(caption = "Tabela 2. Podgląd pierwszych 10 obserwacji") %>% 
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),  full_width = TRUE,font_size = 10)
Tabela 2. Podgląd pierwszych 10 obserwacji
Ref. Limits of Potential Window (V) Lower Limit of Potential Window (V) Upper Limit of Potential Window (V) Potential Window (V) Current Density (A/g) Capacitance (F/g) Specific Surface Area (m^2/g) Charge Transfer Resistance (Rct) (ohm) Equivalent Series Resistance (Rs) (ohm) Electrode Configuration Pore Size (nm) Pore Volume (cm^3/g) Ratio of ID/IG N at% C at% O at% Electrolyte Chemical Formula Electrolyte Ionic Conductivity Electrolyte Concentration (M) Cell Configuration (three/two electrode system)
DOI: 10.1039/c7ta03093b 0 to 0.8 0.0 0.8 0.8 1 680 186.3 NA 7.70 CNF/RGO/moOxNy NA NA 1.450 2.1 NA NA H2SO4 7 1 three-electrode system
DOI: 10.1039/c6ta10933k 0 to 1 0.0 1.0 1.0 1 367 537.0 6.1 1.95 sulfur-doped graphene foam (SGF) NA NA 1.280 0.0 85.6 9.1 KOH 6 6 two-electrode system
DOI: 10.1039/c6ta10933k 0 to 1 0.0 1.0 1.0 2 338 537.0 6.1 1.95 sulfur-doped graphene foam (SGF) NA NA 1.280 0.0 85.6 9.1 KOH 6 6 two-electrode system
DOI: 10.1039/c6ta10933k 0 to 1 0.0 1.0 1.0 5 283 537.0 6.1 1.95 sulfur-doped graphene foam (SGF) NA NA 1.280 0.0 85.6 9.1 KOH 6 6 two-electrode system
DOI: 10.1039/c6ta10933k 0 to 1 0.0 1.0 1.0 10 246 537.0 6.1 1.95 sulfur-doped graphene foam (SGF) NA NA 1.280 0.0 85.6 9.1 KOH 6 6 two-electrode system
DOI: 10.1039/d1se01134k 0 to 0.5 0.0 0.5 0.5 1 872 168.2 NA 0.80 moS2-moO2/G NA NA 1.220 NA NA NA KOH 6 2 three-electrode system
DOI: 10.1039/c8ra08762h -0.4 to 0.2 -0.4 0.2 0.6 1 143 NA NA NA copper anchored boron doped graphene nanosheet (CuBG) 1 NA NA 1.000 NA NA NA H2SO4/PVA NA NA two-electrode system
DOI: 10.1039/c8ra08762h -0.4 to 0.2 -0.4 0.2 0.6 1 306 NA NA NA copper anchored boron doped graphene nanosheet (CuBG) 2 NA NA 1.010 NA NA NA H2SO4/PVA NA NA two-electrode system
DOI: 10.1039/c8ra08762h -0.4 to 0.2 -0.4 0.2 0.6 1 360 NA NA NA copper anchored boron doped graphene nanosheet (CuBG) 3 NA NA 1.033 NA NA NA H2SO4/PVA NA NA two-electrode system
DOI: 10.1039/c8ra08762h -0.4 to 0.2 -0.4 0.2 0.6 1 483 NA NA NA copper anchored boron doped graphene nanosheet (CuBG) 4 NA NA 1.210 NA NA NA H2SO4/PVA NA NA two-electrode system

Zbiór danych zawiera informacje o elektrodach grafitowych testowanych w warunkach laboratoryjnych. Każda obserwacja reprezentuje pojedynczy eksperyment z wyszczególnioną konfiguracją, materiałami i parametrami pomiarowymi.

Przetwarzanie brakujących danych

Analiza przed czyszczeniem

missing_summary <- data %>%
  summarise(across(everything(), ~ sum(is.na(.)))) %>%
  pivot_longer(everything(), names_to = "Zmienna", values_to = "Liczba_NA") %>%  
  mutate(Procent_NA = round(Liczba_NA / nrow(data) * 100, 2)) %>%
  arrange(desc(Procent_NA))    

missing_summary %>%
  kable(caption = "Tabela 3. Procent wartości brakujących przed czyszczeniem",  
        col.names = c("Zmienna", "Liczba NA", "Procent NA (%)")) %>% 
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),  
        full_width = FALSE) %>%
  row_spec(which(missing_summary$Procent_NA > 70), bold = TRUE, color = "white", background = "#d45a57")
Tabela 3. Procent wartości brakujących przed czyszczeniem
Zmienna Liczba NA Procent NA (%)
Charge Transfer Resistance (Rct) (ohm) 786 84.97
Equivalent Series Resistance (Rs) (ohm) 772 83.46
Pore Size (nm) 769 83.14
Pore Volume (cm^3/g) 729 78.81
O at% 703 76.00
C at% 699 75.57
N at% 690 74.59
Ratio of ID/IG 596 64.43
Specific Surface Area (m^2/g) 572 61.84
Electrolyte Ionic Conductivity 99 10.70
Electrolyte Concentration (M) 62 6.70
Electrolyte Chemical Formula 22 2.38
Capacitance (F/g) 17 1.84
Current Density (A/g) 16 1.73
Cell Configuration (three/two electrode system) 14 1.51
Potential Window (V) 5 0.54
Limits of Potential Window (V) 4 0.43
Lower Limit of Potential Window (V) 4 0.43
Upper Limit of Potential Window (V) 4 0.43
Ref. 0 0.00
Electrode Configuration 0 0.00
plot_missing(data, title = "Wizualizacja wartości brakujących")

Analiza wykazała istotne braki w zmiennych: Specific Surface Area (61%), Charge Transfer Resistance (85%), Equivalent Series Resistance (83%), oraz parametrach struktury porowatej i składu chemicznego (69-83% braków).

Charakterystyka typów danych

data.frame(
  Kolumna = names(data),
  Typ = sapply(data, class)
) %>%
  kable(caption = "Tabela 4. Typy danych w kolumnach") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),full_width = FALSE) 
Tabela 4. Typy danych w kolumnach
Kolumna Typ
Ref. Ref. character
Limits of Potential Window (V) Limits of Potential Window (V) character
Lower Limit of Potential Window (V) Lower Limit of Potential Window (V) numeric
Upper Limit of Potential Window (V) Upper Limit of Potential Window (V) numeric
Potential Window (V) Potential Window (V) numeric
Current Density (A/g) Current Density (A/g) numeric
Capacitance (F/g) Capacitance (F/g) numeric
Specific Surface Area (m^2/g) Specific Surface Area (m^2/g) numeric
Charge Transfer Resistance (Rct) (ohm) Charge Transfer Resistance (Rct) (ohm) numeric
Equivalent Series Resistance (Rs) (ohm) Equivalent Series Resistance (Rs) (ohm) numeric
Electrode Configuration Electrode Configuration character
Pore Size (nm) Pore Size (nm) numeric
Pore Volume (cm^3/g) Pore Volume (cm^3/g) numeric
Ratio of ID/IG Ratio of ID/IG numeric
N at% N at% numeric
C at% C at% numeric
O at% O at% numeric
Electrolyte Chemical Formula Electrolyte Chemical Formula character
Electrolyte Ionic Conductivity Electrolyte Ionic Conductivity numeric
Electrolyte Concentration (M) Electrolyte Concentration (M) numeric
Cell Configuration (three/two electrode system) Cell Configuration (three/two electrode system) character
num_vars <- data %>% select(where(is.numeric)) 
data.frame(
  Zmienna = names(num_vars), 
  Brakujące = sapply(num_vars, function(x) sum(is.na(x))),
  Procent_NA = sapply(num_vars, function(x) round(sum(is.na(x))/length(x)*100, 2)) 
) %>%
  kable(caption = "Tabela 5. Zmienne numeryczne - statystyka braków",
        col.names = c("Zmienna", "Liczba NA", "Procent NA (%)")) %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Tabela 5. Zmienne numeryczne - statystyka braków
Zmienna Liczba NA Procent NA (%)
Lower Limit of Potential Window (V) Lower Limit of Potential Window (V) 4 0.43
Upper Limit of Potential Window (V) Upper Limit of Potential Window (V) 4 0.43
Potential Window (V) Potential Window (V) 5 0.54
Current Density (A/g) Current Density (A/g) 16 1.73
Capacitance (F/g) Capacitance (F/g) 17 1.84
Specific Surface Area (m^2/g) Specific Surface Area (m^2/g) 572 61.84
Charge Transfer Resistance (Rct) (ohm) Charge Transfer Resistance (Rct) (ohm) 786 84.97
Equivalent Series Resistance (Rs) (ohm) Equivalent Series Resistance (Rs) (ohm) 772 83.46
Pore Size (nm) Pore Size (nm) 769 83.14
Pore Volume (cm^3/g) Pore Volume (cm^3/g) 729 78.81
Ratio of ID/IG Ratio of ID/IG 596 64.43
N at% N at% 690 74.59
C at% C at% 699 75.57
O at% O at% 703 76.00
Electrolyte Ionic Conductivity Electrolyte Ionic Conductivity 99 10.70
Electrolyte Concentration (M) Electrolyte Concentration (M) 62 6.70
cat_vars <- data %>% select(where(~ is.character(.) | is.factor(.)))
data.frame(
  Zmienna = names(cat_vars),   
  Brakujące = sapply(cat_vars, function(x) sum(is.na(x))), 
  Kategorie = sapply(cat_vars, function(x) n_distinct(x, na.rm = TRUE)) 
) %>%
  kable(caption = "Tabela 6. Zmienne kategoryczne - statystyka",
        col.names = c("Zmienna", "Liczba NA", "Liczba kategorii")) %>% 
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))   
Tabela 6. Zmienne kategoryczne - statystyka
Zmienna Liczba NA Liczba kategorii
Ref. Ref. 0 198
Limits of Potential Window (V) Limits of Potential Window (V) 4 63
Electrode Configuration Electrode Configuration 0 353
Electrolyte Chemical Formula Electrolyte Chemical Formula 22 23
Cell Configuration (three/two electrode system) Cell Configuration (three/two electrode system) 14 2

Czyszczenia danych

Etap 1: Usunięcie kolumn z wysokim odsetkiem braków

Kolumny zawierające ponad 70% wartości brakujących zostały usunięte z analizy.

  • Nadmierna imputacja zniekształciłaby rozkłady i zafałszowala wyniki analiz statystycznych,
  • Usunięcie wierszy zmniejszyłoby zbiór do poniżej 200 obserwacji.

Etap 2: Wyliczenie wartości dla okna potencjału

Wykorzystano relacje między dolną granicą, górną granicą i szerokością okna potencjalu.

Etap 3: Imputacja pozostałych braków

  • Zmienne numeryczne: mediana (odporna na wartości odstające),
  • Zmienne kategoryczne: moda - najczęstsza kategoria.

Etap 4: Usunięcie obserwacji bez zmiennej docelowej

Wiersze bez wartości Capacitance zostały usunięte.

Implementacja czyszczenia

na_before <- sapply(data, function(x) sum(is.na(x)))

data_fix <- data %>%
  mutate(
    `Potential Window (V)` = coalesce(`Potential Window (V)`,
                                      `Upper Limit of Potential Window (V)` - `Lower Limit of Potential Window (V)`),
    `Upper Limit of Potential Window (V)` = coalesce(`Upper Limit of Potential Window (V)`,
                                                    `Lower Limit of Potential Window (V)` + `Potential Window (V)`),
    `Lower Limit of Potential Window (V)` = coalesce(`Lower Limit of Potential Window (V)`,
                                                    `Upper Limit of Potential Window (V)` - `Potential Window (V)`),
    `Limits of Potential Window (V)` = coalesce(`Limits of Potential Window (V)`,
                                                paste(`Lower Limit of Potential Window (V)`, "to", `Upper Limit of Potential Window (V)`)))

na_after <- sapply(data_fix, function(x) sum(is.na(x)))

cols_window <- c("Potential Window (V)", "Upper Limit of Potential Window (V)", 
                 "Lower Limit of Potential Window (V)", "Limits of Potential Window (V)")
filled_summary <- data.frame(
  Kolumna = cols_window,
  NA_przed = na_before[cols_window],
  NA_po = na_after[cols_window],
  Uzupełnione = na_before[cols_window] - na_after[cols_window]
)

filled_summary %>%
  kable(caption = "Tabela 7. Uzupełnienie wartości okna potencjału") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Tabela 7. Uzupełnienie wartości okna potencjału
Kolumna NA_przed NA_po Uzupełnione
Potential Window (V) Potential Window (V) 5 4 1
Upper Limit of Potential Window (V) Upper Limit of Potential Window (V) 4 4 0
Lower Limit of Potential Window (V) Lower Limit of Potential Window (V) 4 4 0
Limits of Potential Window (V) Limits of Potential Window (V) 4 0 4
data_clean <- data_fix %>%
  select(-`Ref.`, -`Electrode Configuration`) %>%
  select(where(~ mean(is.na(.)) < 0.7))

rows_before <- nrow(data_clean)
data_clean <- data_clean %>% filter(!is.na(`Capacitance (F/g)`))
rows_after <- nrow(data_clean)
removed_rows <- rows_before - rows_after

cat("\nUsunięto wierszy bez zmiennej docelowej:", removed_rows, "\n")
## 
## Usunięto wierszy bez zmiennej docelowej: 17
cat("Pozostało obserwacji:", rows_after, "\n\n")
## Pozostało obserwacji: 908
data_imputed <- data_clean

num_cols <- names(data_imputed)[sapply(data_imputed, is.numeric)]
for (col in num_cols) {
  data_imputed[[col]][is.na(data_imputed[[col]])] <- median(data_imputed[[col]], na.rm = TRUE)
}

cat_cols <- names(data_imputed)[sapply(data_imputed, function(x) is.character(x) | is.factor(x))]
for (col in cat_cols) {
  mode_val <- names(sort(table(data_imputed[[col]]), decreasing = TRUE))[1]
  data_imputed[[col]][is.na(data_imputed[[col]])] <- mode_val
}

rename_map <- c(
  "Specific Surface Area (m^2/g)" = "SSA",
  "Capacitance (F/g)" = "Capacitance",
  "Cell Configuration (three/two electrode system)" = "Cell_Config"
)

existing_names <- names(rename_map)[names(rename_map) %in% colnames(data_imputed)]
data_imputed <- data_imputed %>%
  rename(!!!setNames(existing_names, rename_map[existing_names]))

cat("Finalne wymiary zbioru po czyszczeniu:", nrow(data_imputed), "×", ncol(data_imputed), "\n")
## Finalne wymiary zbioru po czyszczeniu: 908 × 12

Rozmiar zbioru i podstawowe statystyki

Podsumowanie zbioru po czyszczeniu

plot_missing(data_imputed, title = "Kompletność danych po czyszczeniu")

cat("\nLiczba obserwacji:", nrow(data_imputed), "\n")
## 
## Liczba obserwacji: 908
cat("Liczba zmiennych:", ncol(data_imputed), "\n\n")
## Liczba zmiennych: 12
num_vars <- data_imputed %>% select(where(is.numeric))
data.frame(
  Zmienna = names(num_vars),
  Brakujące = sapply(num_vars, function(x) sum(is.na(x)))
) %>%
  kable(caption = "Tabela 8. Zmienne numeryczne po czyszczeniu",
        col.names = c("Zmienna", "Liczba NA")) %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Tabela 8. Zmienne numeryczne po czyszczeniu
Zmienna Liczba NA
Lower Limit of Potential Window (V) Lower Limit of Potential Window (V) 0
Upper Limit of Potential Window (V) Upper Limit of Potential Window (V) 0
Potential Window (V) Potential Window (V) 0
Current Density (A/g) Current Density (A/g) 0
Capacitance Capacitance 0
SSA SSA 0
Ratio of ID/IG Ratio of ID/IG 0
Electrolyte Ionic Conductivity Electrolyte Ionic Conductivity 0
Electrolyte Concentration (M) Electrolyte Concentration (M) 0
cat_vars <- data_imputed %>% select(where(~ is.character(.) | is.factor(.)))
data.frame(
  Zmienna = names(cat_vars),
  Brakujące = sapply(cat_vars, function(x) sum(is.na(x))),
  Kategorie = sapply(cat_vars, function(x) n_distinct(x, na.rm = TRUE))
) %>%
  kable(caption = "Tabela 9. Zmienne kategoryczne po czyszczeniu",
        col.names = c("Zmienna", "Liczba NA", "Liczba kategorii")) %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Tabela 9. Zmienne kategoryczne po czyszczeniu
Zmienna Liczba NA Liczba kategorii
Limits of Potential Window (V) Limits of Potential Window (V) 0 60
Electrolyte Chemical Formula Electrolyte Chemical Formula 0 22
Cell_Config Cell_Config 0 2

Statystyki opisowe

num_data <- data_imputed %>% select(where(is.numeric))
num_data %>%
  psych::describe() %>%
  select(n, mean, sd, median, min, max, range, skew, kurtosis) %>%
  kable(digits = 2, caption = "Tabela 10. Statystyki opisowe zmiennych numerycznych") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = TRUE,
                font_size = 10)
Tabela 10. Statystyki opisowe zmiennych numerycznych
n mean sd median min max range skew kurtosis
Lower Limit of Potential Window (V) 908 -0.23 0.37 0.00 -1.10 0.20 1.30 -1.33 0.15
Upper Limit of Potential Window (V) 908 0.63 0.45 0.60 -0.20 3.50 3.70 1.57 7.32
Potential Window (V) 908 0.86 0.35 0.85 0.40 3.50 3.10 2.70 14.64
Current Density (A/g) 908 5.86 13.36 2.00 0.05 200.00 199.95 8.68 106.15
Capacitance 908 415.50 447.53 260.25 1.40 3344.08 3342.68 2.75 9.46
SSA 908 257.81 360.43 159.97 8.90 2400.00 2391.10 3.64 13.97
Ratio of ID/IG 908 1.08 0.26 1.05 0.12 2.90 2.78 3.16 20.53
Electrolyte Ionic Conductivity 908 5.83 1.31 6.00 1.00 8.00 7.00 -1.67 2.28
Electrolyte Concentration (M) 908 2.49 2.16 1.00 0.10 6.00 5.90 0.89 -1.02

Finalny zbiór danych zawiera 908 obserwacji i 12 zmiennych, w tym 9 numerycznych i 3 kategoryczne, wszystkie bez braków. Zmienne numeryczne charakteryzują się dużą zmiennością i asymetrią rozkładów - wiele z nich jest prawoskośnych, co oznacza, że większość wartości skupia się w dolnym zakresie, a pojedyncze obserwacje mają znacznie większe wartości. Podnosi to średnią powyżej mediany. Na przykład pojemność elektrod waha się od 1.4 do 3344 F/g, ze średnią 415.5 i medianą 260.3 F/g, co wskazuje na wyraźną prawostronną asymetrię.

Natężenie prądu w stosunku do masy (Current density) ma ekstremalną wartość kurtozy - 106.15. Kurtoza opisuje „spiczastość” rozkładu i koncentrację obserwacji wokół średniej. Ta wysoka wartość ozncza, że parametr ma wiele ekstremalnych punktów, które mogą istotnie wpływać na analizę i modele ML. Potencjał elektrody wykazuje bardziej symetryczne rozkłady, a przewodność elektrolitu i jego stężenie są względnie stabilne. Prewodniść elektrolit est zmienną liczbową przyjmującą ogrnaniconą liczbę wartości całkowitych (od 1 do 8). W przypadku stężenia elektrolitu obserwuje się największą zmienność, ze standardowym odchyleniem wynoszącym 2. Wnioski dla analizy ML są takie, że należy uwzględnić nieliniowości, stosować modele odporne na wartości odstające oraz ewentualnie transformacje danych (np. log), szczególnie dla cech takich jak Capacitance czy Current Density.

Szczegółowa analiza wartości atrybutów

Rozkłady zmiennych numerycznych

num_cols <- data_imputed %>% select(where(is.numeric)) %>% names()

plots_list <- lapply(num_cols, function(col) {
  ggplot(data_imputed, aes(x = .data[[col]])) +
    geom_histogram(bins = 30, fill = "steelblue", color = "white", alpha = 0.8) +
    geom_density(aes(y = after_stat(count)), color = "red", linewidth = 1) +
    labs(title = col, x = "", y = "Częstość") +
    theme_minimal() +
    theme(plot.title = element_text(size = 10, face = "bold"))
})

gridExtra::grid.arrange(grobs = plots_list, ncol = 3, 
                       top = "Rozkłady zmiennych numerycznych")

Rozkłady zmiennych kategorycznych

cat_cols <- data_imputed %>% select(where(~ is.character(.) | is.factor(.))) %>% names()

for (col in cat_cols) {
  p <- ggplot(data_imputed, aes(x = .data[[col]])) + 
    geom_bar(fill = "steelblue", alpha = 0.8) +
    labs(title = paste("Rozkład:", col), x = "", y = "Liczba obserwacji") + 
    theme_minimal() + 
    theme(axis.text.x = element_text(angle = 45, hjust = 1, size = 8))
  print(p)
}

Zmienna opisująca limity przedziałów wartości potencjału jest zmienną kategoryczną o bardzo dużej liczbie kategorii i w dużej mierze powiela informacje zawarte już w zmiennych liczbowych. Nie będzie uwzględniana w dalszej analizie ani w modelach ML.

Pozostałe dwie zmienne kategoryczne mają akceptowalną liczbę kategorii - w danych występują 22 różne typy elektrolitów, z wyraźną dominacją KOH oraz H₂SO₄, oraz dwa typy konfiguracji konduktora. Te zmienne zachowujemy do dalszej analizy.

Test normalności rozkładów

check_normality <- function(data) {
  num_cols <- names(data %>% select(where(is.numeric)))
  results <- list()
  for (col in num_cols) {
    if (length(na.omit(data[[col]])) > 3) {
      p_val <- shapiro.test(data[[col]])$p.value
      results[[col]] <- ifelse(p_val > 0.05, "Rozkład Normalny", "Rozkład niezgodny z rozkładem normalnym")
    }
  }
  return(results)   
}  

normality_results <- check_normality(data_imputed)

data.frame(
  Zmienna = names(normality_results),
  Test_Shapiro_Wilka = unlist(normality_results)  
) %>%
  kable(caption = "Tabela 11. Wyniki testu normalności Shapiro-Wilka",
        col.names = c("Zmienna", "Wynik testu (α = 0.05)")) %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))  
Tabela 11. Wyniki testu normalności Shapiro-Wilka
Zmienna Wynik testu (α = 0.05)
Lower Limit of Potential Window (V) Lower Limit of Potential Window (V) Rozkład niezgodny z rozkładem normalnym
Upper Limit of Potential Window (V) Upper Limit of Potential Window (V) Rozkład niezgodny z rozkładem normalnym
Potential Window (V) Potential Window (V) Rozkład niezgodny z rozkładem normalnym
Current Density (A/g) Current Density (A/g) Rozkład niezgodny z rozkładem normalnym
Capacitance Capacitance Rozkład niezgodny z rozkładem normalnym
SSA SSA Rozkład niezgodny z rozkładem normalnym
Ratio of ID/IG Ratio of ID/IG Rozkład niezgodny z rozkładem normalnym
Electrolyte Ionic Conductivity Electrolyte Ionic Conductivity Rozkład niezgodny z rozkładem normalnym
Electrolyte Concentration (M) Electrolyte Concentration (M) Rozkład niezgodny z rozkładem normalnym

Wszystkie zmienne numeryczne wykazują istotne odchylenia od rozkładu normalnego (test Shapiro-Wilka, p < 0.05). Rozkłady charakteryzują się silną skośnością i wysoką kurtozą. W konsekwencji w dalszej analizie zastosowano współczynnik korelacji Spearmana.

Korelacje między zmiennymi

Korelacje zmiennych numerycznych

cor_matrix_spear <- cor(num_data, use = "pairwise.complete.obs", method = "spearman")

corrplot(cor_matrix_spear, method = "color", type = "upper", 
         tl.cex = 0.8, tl.col = "black", addCoef.col = "black", 
         number.cex = 0.6,
         title = "Macierz korelacji Spearmana",
         mar = c(0,0,2,0))

cor_spear_df <- as.data.frame(cor_matrix_spear)
cor_spear_df$Var1 <- rownames(cor_spear_df)
cor_long <- cor_spear_df %>%
  pivot_longer(-Var1, names_to = "Var2", values_to = "Corr") %>%
  filter(Var1 != Var2) %>%
  mutate(AbsCorr = abs(Corr)) %>%
  arrange(desc(AbsCorr)) %>%
  distinct(AbsCorr, .keep_all = TRUE) %>%
  head(10)

cor_long %>%
  select(Var1, Var2, Corr) %>%
  kable(digits = 3, caption = "Tabela 12. Najsilniejsze korelacje Spearmana",
        col.names = c("Zmienna 1", "Zmienna 2", "Współczynnik ρ")) %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Tabela 12. Najsilniejsze korelacje Spearmana
Zmienna 1 Zmienna 2 Współczynnik ρ
Lower Limit of Potential Window (V) Upper Limit of Potential Window (V) 0.620
Potential Window (V) Capacitance -0.437
Upper Limit of Potential Window (V) Electrolyte Concentration (M) -0.424
Upper Limit of Potential Window (V) Potential Window (V) 0.359
Lower Limit of Potential Window (V) Potential Window (V) -0.334
Upper Limit of Potential Window (V) Electrolyte Ionic Conductivity 0.301
Lower Limit of Potential Window (V) Electrolyte Ionic Conductivity 0.237
Potential Window (V) SSA 0.189
SSA Electrolyte Concentration (M) 0.184
Current Density (A/g) Electrolyte Concentration (M) 0.183

Analiza korelacji wykazała słabe zależności między zmiennymi. Najsilniejsze korelacje obserwowano między parametrami dotyczącymi przedzialów potencjału, ponieważ są one od siebie nauralnie zależne. Pojemność elektrody (Capacitance) wykazuje najsilniejszą korelację z przedziałem potencjału (-0.44). Ujemny znak korelacji wskazuje, że szersze przedziały potencjałusą związane z niższą pojemnością w porównaniu do wąskich przedziałów.

Korelacje zmiennych kategorycznych

cat_data <- data_imputed %>% select(where(~ is.character(.) | is.factor(.)))

if (ncol(cat_data) > 1) {
  cramer_matrix <- matrix(NA, ncol = ncol(cat_data), nrow = ncol(cat_data))
  colnames(cramer_matrix) <- rownames(cramer_matrix) <- names(cat_data)
  
  for (i in 1:ncol(cat_data)) {
    for (j in 1:ncol(cat_data)) {
      if (i != j) {
        cramer_matrix[i, j] <- DescTools::CramerV(table(cat_data[[i]], cat_data[[j]]))
      } else {
        cramer_matrix[i, j] <- 1
      }
    }
  }
  
  corrplot(cramer_matrix, method = "color", type = "upper", 
           tl.cex = 0.8, tl.col = "black", addCoef.col = "black", 
           number.cex = 0.7,
           title = "Macierz asocjacji Cramer's V",
           mar = c(0,0,2,0))
}

Zmienne kategoryczne porównane zostały z użyciem Cramer’s V - miara siły związku między zmiennymi kategorycznymi. Zmienne kategoryczne opisują materiały użyte do ekspermentu oraz jego konfigurację ogniwa - wykazują więc dość duże korelacje między sobą.

Zależności kategoryczne-numeryczne

Aby zbadać wpływ zmiennych kategorycznych (typ elektrolitu, konfiguracja ogniwa) na zmienne numeryczne, standardowo stosuje się ANOVA i Eta-squared.ANOVA wymaga spełnienia założeń:

Normalność reszt w każdej grupie (test Shapiro-Wilka) Homogeniczność wariancji między grupami (test Levene’a)

Sprawdzenie:

cat_cols <- names(data_imputed %>% select(where(~ is.character(.) | is.factor(.))))
num_cols <- names(data_imputed %>% select(where(is.numeric)))

cat("### Test założeń ANOVA dla wybranych par zmiennych:\n\n")
## ### Test założeń ANOVA dla wybranych par zmiennych:
test_pairs <- list(
  c("Electrolyte Chemical Formula", "Capacitance"),
  c("Cell_Config", "Capacitance"),
  c("Electrolyte Chemical Formula", "Current Density (A/g)"),
  c("Cell_Config", "SSA"),
  c("Electrolyte Chemical Formula", "Potential Window (V)")
)

violations <- 0
total <- 0

for (pair in test_pairs) {
  cat_var <- pair[1]
  num_var <- pair[2]
  
  if (cat_var %in% cat_cols && num_var %in% num_cols) {
    total <- total + 1
    
    #ANOVA
    anova_model <- aov(data_imputed[[num_var]] ~ data_imputed[[cat_var]])
    
    #Test normalności reszt
    residuals_vals <- residuals(anova_model)
    if (length(residuals_vals) > 5000) {
      residuals_vals <- sample(residuals_vals, 1000)
    }
    shapiro_p <- shapiro.test(residuals_vals)$p.value
    
    #Test homogeniczności wariancji
    levene_test <- car::leveneTest(data_imputed[[num_var]] ~ as.factor(data_imputed[[cat_var]]))
    levene_p <- levene_test$`Pr(>F)`[1]
    
    #Sprawdzenie założeń
    norm_ok <- shapiro_p > 0.05
    homo_ok <- levene_p > 0.05
    
    cat("**", pair[1], "→", pair[2], ":**\n")
    cat("- Normalność:", ifelse(norm_ok, "Tak", "Nie"), "(p =", round(shapiro_p, 4), ")\n")
    cat("- Homogeniczność:", ifelse(homo_ok, "Tak", "Nie"), "(p =", round(levene_p, 4), ")\n\n")
    
    if (!norm_ok || !homo_ok) {
      violations <- violations + 1
    }
  }
}
## ** Electrolyte Chemical Formula → Capacitance :**
## - Normalność: Nie (p = 0 )
## - Homogeniczność: Nie (p = 0 )
## 
## ** Cell_Config → Capacitance :**
## - Normalność: Nie (p = 0 )
## - Homogeniczność: Nie (p = 0 )
## 
## ** Electrolyte Chemical Formula → Current Density (A/g) :**
## - Normalność: Nie (p = 0 )
## - Homogeniczność: Nie (p = 0 )
## 
## ** Cell_Config → SSA :**
## - Normalność: Nie (p = 0 )
## - Homogeniczność: Nie (p = 0 )
## 
## ** Electrolyte Chemical Formula → Potential Window (V) :**
## - Normalność: Nie (p = 0 )
## - Homogeniczność: Nie (p = 0 )
cat("**", violations, "z", total, "par (", round(violations/total*100, 1), "%) Nie spełnia założeń ANOVA.**\n\n")
## ** 5 z 5 par ( 100 %) Nie spełnia założeń ANOVA.**

Interaktywne wykresy

Wykres punktowy z dropdown

num_vars <- data_imputed %>% select(where(is.numeric)) %>% names()
num_vars <- num_vars[num_vars != "Capacitance"]

if (length(num_vars) > 0) {
  tooltip_list <- lapply(num_vars, function(var) {
    do.call(paste, c(lapply(colnames(data_imputed), function(col) {
      paste0(col, ": ", data_imputed[[col]])
    }), sep = "<br>"))
  })
  names(tooltip_list) <- num_vars
  
  axis_ranges <- lapply(num_vars, function(var) {
    x_vals <- data_imputed[[var]]
    x_min <- min(x_vals, na.rm = TRUE)
    x_max <- max(x_vals, na.rm = TRUE)
    x_range <- x_max - x_min
    margin <- x_range * 0.1
    list(min = x_min - margin, max = x_max + margin)
  })
  names(axis_ranges) <- num_vars
  
  create_scatter <- function(x_var) {
    p <- ggplot(data_imputed, aes(x = .data[[x_var]], y = Capacitance, text = tooltip_list[[x_var]])) +
      geom_point(alpha = 0.6, size = 2, color = "#2c7bb6") +
      geom_smooth(method = "loess", se = TRUE, color = "red", linewidth = 1) +
      labs(title = "Zależność pojemności od wybranej zmiennej", 
           x = x_var, 
           y = "Capacitance (F/g)") +
      theme_minimal()
    
    ggplotly(p, tooltip = "text")
  }
  
  scatter_plot <- create_scatter(num_vars[1])
  
  scatter_plot <- scatter_plot %>%
    layout(
      updatemenus = list(
        list(
          type = "dropdown",
          buttons = lapply(num_vars, function(var) {
            list(
              method = "update",
              args = list(
                list(x = list(data_imputed[[var]]), text = list(tooltip_list[[var]])),
                list(
                  xaxis = list(
                    title = var,
                    range = c(axis_ranges[[var]]$min, axis_ranges[[var]]$max)
                  )
                )
              ),
              label = var
            )
          }),
          direction = "down",
          showactive = TRUE,
          x = 0.1,
          xanchor = "left",
          y = 1.15,
          yanchor = "top"
        )
      )
    )
  
  scatter_plot
}

Interaktywny wykres punktowy umożliwia eksplorację zależności między pojemnością a pozostałymi zmiennymi numerycznymi. Funkcja hover wyświetla wszystkie parametry danej elektrody po najechaniu na punkt, co ułatwia identyfikację właściwości materaiłów. Wykresy potwierdzają wysokie rozproszenie danych i brak prostych liniowych zależności.

Przygotowanie danych do modelowania ML

Analiza wartości odstających

Przed budową modeli należy zbadać wartości odstające (outliers), które mogą obniżać dokładność modelum który uczy się anomalii zamiast prawdziwych wzorców.

W analizie brana będzie pod uwagę zmienna Capacitance (zmienna docelowa): - ma dużw wartości kurtozy oraz skośności - zakres 1.4-3344 F/g, a mediana to 260 F/g, więc pojedyncze wartości >1000 F/g mogą zdominować model

oraz Current Density: - ma bardzo wysoką kurtozę (106.15) - mediana to 2 A/g, ale wartości odstajace mają >50 A/g

#statystyki IQR
num_cols <- names(data_imputed %>% select(where(is.numeric)))

outlier_stats <- data.frame(
  Zmienna = character(),
  Q1 = numeric(),
  Q3 = numeric(),
  IQR = numeric(),
  Dolny_prog = numeric(),
  Gorny_prog = numeric(),
  Liczba_outliers = integer(),
  Procent_outliers = numeric(),
  stringsAsFactors = FALSE
)

for (col in num_cols) {
  q1 <- quantile(data_imputed[[col]], 0.25, na.rm = TRUE)
  q3 <- quantile(data_imputed[[col]], 0.75, na.rm = TRUE)
  iqr <- q3 - q1
  lower <- q1 - 1.5 * iqr
  upper <- q3 + 1.5 * iqr
  
  outliers <- data_imputed[[col]] < lower | data_imputed[[col]] > upper
  n_outliers <- sum(outliers, na.rm = TRUE)
  pct_outliers <- round(n_outliers / nrow(data_imputed) * 100, 2)
  
  outlier_stats <- rbind(outlier_stats, data.frame(
    Zmienna = col,
    Q1 = q1,
    Q3 = q3,
    IQR = iqr,
    Dolny_prog = lower,
    Gorny_prog = upper,
    Liczba_outliers = n_outliers,
    Procent_outliers = pct_outliers
  ))
}

outlier_stats %>%
  arrange(desc(Procent_outliers)) %>%
  kable(digits = 2, caption = "Tabela 17. Statystyki wartości odstających (metoda IQR 1.5×)") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed")) %>%
  row_spec(which(outlier_stats$Procent_outliers > 5), bold = TRUE, color = "white", background = "#e74c3c")
Tabela 17. Statystyki wartości odstających (metoda IQR 1.5×)
Zmienna Q1 Q3 IQR Dolny_prog Gorny_prog Liczba_outliers Procent_outliers
25%5 SSA 159.97 159.97 0.00 159.97 159.97 340 37.44
25%6 Ratio of ID/IG 1.05 1.05 0.00 1.05 1.05 316 34.80
25% Lower Limit of Potential Window (V) -0.30 0.00 0.30 -0.75 0.45 156 17.18
25%7 Electrolyte Ionic Conductivity 6.00 7.00 1.00 4.50 8.50 121 13.33
25%3 Current Density (A/g) 1.00 5.00 4.00 -5.00 11.00 96 10.57
25%4 Capacitance 148.60 509.85 361.25 -393.28 1051.72 67 7.38
25%1 Upper Limit of Potential Window (V) 0.40 0.80 0.40 -0.20 1.40 23 2.53
25%2 Potential Window (V) 0.60 1.00 0.40 0.00 1.60 19 2.09
25%8 Electrolyte Concentration (M) 1.00 6.00 5.00 -6.50 13.50 0 0.00

Zmienne, w których udział wartości odstających przekracza 5%, wymagają dodatkowej analizy, ponieważ mogą zaburzać wyniki modelowania. W metodzie IQR za obserwacje odstające uznaje się te, które leżą poniżej Q1 − 1.5xIQR lub powyżej Q3 + 1.5xQR, co pozwala wyznaczyć wartości nietypowe .

Usuwanie ekstremalnych wartości odstających

Zastosujemy proste obcięcie wartości ekstremalnych - usunięcie górne oraz dolne 0.5% wartości dla wszystkich zmiennych z >10% outlierów: SSA, Ratio of ID/IG, Lower Limit of Potential Window (V), Electrolyte Ionic Conductivity, Current Density (A/g).

#Zmienne do obcięcia (>10% outlierów)
vars_to_clip <- c("SSA", "Ratio of ID/IG", "Lower Limit of Potential Window (V)", 
                  "Electrolyte Ionic Conductivity", "Current Density (A/g)")

data_cleaned <- data_imputed

cat("### Progi obcięcia (0.5% - 99.5%):\n")
## ### Progi obcięcia (0.5% - 99.5%):
for (var in vars_to_clip) {
  lower <- quantile(data_imputed[[var]], 0.005, na.rm = TRUE)
  upper <- quantile(data_imputed[[var]], 0.995, na.rm = TRUE)
  data_cleaned[[var]] <- pmin(pmax(data_cleaned[[var]], lower), upper)
  cat(var, ": [", round(lower, 2), ",", round(upper, 2), "]\n")
}
## SSA : [ 8.9 , 2211.1 ]
## Ratio of ID/IG : [ 0.16 , 2.81 ]
## Lower Limit of Potential Window (V) : [ -1.05 , 0.2 ]
## Electrolyte Ionic Conductivity : [ 1.54 , 8 ]
## Current Density (A/g) : [ 0.1 , 89.3 ]

Finalne przygotowanie danych

model_data <- data_cleaned %>%
  drop_na() %>%
  select(-`Limits of Potential Window (V)`) %>%  #Usunięcie zmiennej tekstowej mającje już pokreyce w zmiennych numerycznych
  mutate(
    `Electrolyte Chemical Formula` = as.factor(`Electrolyte Chemical Formula`),
    Cell_Config = as.factor(Cell_Config)
  )

Podział na zbiór treningowy i testowy

set.seed(123)
train_index <- createDataPartition(model_data$Capacitance, p = 0.8, list = FALSE)
train_data <- model_data[train_index, ]
test_data <- model_data[-train_index, ]

X_train <- train_data %>% select(-Capacitance)
y_train <- train_data$Capacitance
X_test <- test_data %>% select(-Capacitance)
y_test <- test_data$Capacitance

cat("Rozmiar zbioru treningowego:", nrow(train_data), "obserwacji\n")
## Rozmiar zbioru treningowego: 728 obserwacji
cat("Rozmiar zbioru testowego:", nrow(test_data), "obserwacji\n")
## Rozmiar zbioru testowego: 180 obserwacji

Dane podzielono stratyfikowanym losowaniem na zbiór treningowy (80%) i testowy (20%). Stratyfikacja zapewnia podobne rozkłady zmiennej docelowej w obu zbiorach.

Predykcja pojemności i wyjaśnialna sztuczna inteligencja

Model 1: Random Forest

ctrl <- trainControl(method = "cv", number = 5)
grid <- expand.grid(mtry = c(2, 4, 6))

rf_model <- train(
  x = X_train,
  y = y_train,
  method = "rf",
  tuneGrid = grid,
  trControl = ctrl,
  ntree = 100,
  importance = TRUE
)

cat("Najlepszy model Random Forest:\n")
## Najlepszy model Random Forest:
print(rf_model$bestTune)
##   mtry
## 2    4
rf_predictions <- predict(rf_model, newdata = X_test)
rf_metrics <- postResample(pred = rf_predictions, obs = y_test)

cat("\nMetryki na zbiorze testowym:\n")
## 
## Metryki na zbiorze testowym:
print(rf_metrics)
##        RMSE    Rsquared         MAE 
## 237.8766768   0.7007953 137.3859221
rf_results_df <- data.frame(
  Observed = y_test,
  Predicted = rf_predictions
)

ggplot(rf_results_df, aes(x = Observed, y = Predicted)) +
  geom_point(alpha = 0.6, color = "#2c7bb6", size = 2) +
  geom_abline(intercept = 0, slope = 1, linetype = "dashed", color = "red", linewidth = 1) +
  labs(title = "Predykcje Random Forest",
       x = "Rzeczywista pojemność (F/g)",
       y = "Przewidywana pojemność (F/g)") +
  theme_minimal()

varImp(rf_model)$importance %>%
  as.data.frame() %>%
  tibble::rownames_to_column("Feature") %>%
  arrange(desc(Overall)) %>%
  head(10) %>%
  kable(digits = 2, caption = "Tabela 14. Ważność cech w modelu Random Forest",
        col.names = c("Cecha", "Ważność")) %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Tabela 14. Ważność cech w modelu Random Forest
Cecha Ważność
Potential Window (V) 100.00
Ratio of ID/IG 91.71
Electrolyte Concentration (M) 71.03
Upper Limit of Potential Window (V) 52.63
Electrolyte Chemical Formula 49.34
SSA 27.70
Lower Limit of Potential Window (V) 24.71
Electrolyte Ionic Conductivity 9.57
Current Density (A/g) 6.17
Cell_Config 0.00

Model Random Forest z optymalnym parametrem mtry=4 osiągnął R² = 0.7 i RMSE = 237 F/g na zbiorze testowym. Wyniki wskazują na umiarkowaną zdolność predykcyjną. Najważniejszymi cechami są przedział wartości potencjału i natężenie podzielone przez masę elektrody (Current density).

Model 2: Support Vector Regression

all_X <- rbind(X_train, X_test)
dummyVars_model <- dummyVars(~ ., data = all_X, fullRank = TRUE)
all_X_dummy <- predict(dummyVars_model, newdata = all_X)

train_rows <- 1:nrow(X_train)
X_train_dummy <- all_X_dummy[train_rows, ]
X_test_dummy <- all_X_dummy[-train_rows, ]

svr_grid <- expand.grid(
  sigma = c(0.01, 0.1, 1),
  C = c(1, 10, 100)
)

svr_control <- trainControl(method = "cv", number = 3, verboseIter = FALSE)

set.seed(123)
svr_model <- train(
  x = X_train_dummy,
  y = y_train,
  method = "svmRadial",
  trControl = svr_control,
  tuneGrid = svr_grid,
  preProcess = c("center", "scale")
)

svr_predictions <- predict(svr_model, newdata = X_test_dummy)
svr_metrics <- postResample(pred = svr_predictions, obs = y_test)

cat("**Najlepsze parametry SVR:** sigma =", svr_model$bestTune$sigma, ", C =", svr_model$bestTune$C, "\n\n")
## **Najlepsze parametry SVR:** sigma = 1 , C = 100
cat("**Metryki na zbiorze testowym:**\n")
## **Metryki na zbiorze testowym:**
cat("- RMSE:", round(svr_metrics["RMSE"], 2), "\n")
## - RMSE: 376.17
cat("- R²:", round(svr_metrics["Rsquared"], 3), "\n")
## - R²: 0.314
cat("- MAE:", round(svr_metrics["MAE"], 2), "\n")
## - MAE: 185.02

Model 3: Elastic Net

elastic_grid <- expand.grid(
  alpha = c(0, 0.5, 1),
  lambda = 10^seq(-3, 1, length = 10)
)

elastic_control <- trainControl(method = "cv", number = 5, verboseIter = FALSE)

set.seed(123)
elastic_model <- train(
  x = X_train,
  y = y_train,
  method = "glmnet",
  trControl = elastic_control,
  tuneGrid = elastic_grid,
  preProcess = c("center", "scale")
)

elastic_predictions <- predict(elastic_model, newdata = X_test)
elastic_metrics <- postResample(pred = elastic_predictions, obs = y_test)

cat("**Najlepsze parametry Elastic Net:** alpha =", elastic_model$bestTune$alpha, ", lambda =", elastic_model$bestTune$lambda, "\n\n")
## **Najlepsze parametry Elastic Net:** alpha = 0 , lambda = 10
cat("**Metryki na zbiorze testowym:**\n")
## **Metryki na zbiorze testowym:**
cat("- RMSE:", round(elastic_metrics["RMSE"], 2), "\n")
## - RMSE: 395.66
cat("- R²:", round(elastic_metrics["Rsquared"], 3), "\n")
## - R²: 0.169
cat("- MAE:", round(elastic_metrics["MAE"], 2), "\n")
## - MAE: 263.68

Porównanie modeli

comparison_df <- data.frame(
  Model = c("Random Forest", "SVR", "Elastic Net"),
  RMSE = c(rf_metrics["RMSE"], svr_metrics["RMSE"], elastic_metrics["RMSE"]),
  R2 = c(rf_metrics["Rsquared"], svr_metrics["Rsquared"], elastic_metrics["Rsquared"]),
  MAE = c(rf_metrics["MAE"], svr_metrics["MAE"], elastic_metrics["MAE"])
)

comparison_df %>%
  kable(digits = 2, caption = "Tabela 15. Porównanie modeli predykcyjnych",
        col.names = c("Model", "RMSE", "R²", "MAE")) %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Tabela 15. Porównanie modeli predykcyjnych
Model RMSE MAE
Random Forest 237.88 0.70 137.39
SVR 376.17 0.31 185.02
Elastic Net 395.66 0.17 263.68
comparison_plot_df <- data.frame(
  Observed = rep(y_test, 3),
  Predicted = c(rf_predictions, svr_predictions, elastic_predictions),
  Model = rep(c("Random Forest", "SVR", "Elastic Net"), each = length(y_test))
)

ggplot(comparison_plot_df, aes(x = Observed, y = Predicted, color = Model)) +
  geom_point(alpha = 0.5, size = 2) +
  geom_abline(intercept = 0, slope = 1, linetype = "dashed", color = "black") +
  facet_wrap(~Model) +
  labs(title = "Porównanie predykcji trzech modeli",
       x = "Rzeczywista pojemność (F/g)",
       y = "Przewidywana pojemność (F/g)") +
  theme_minimal() +
  theme(legend.position = "bottom")

Random Forest okazał się najlepszym modelem spośród trzech testowanych algorytmów, osiągając najniższy RMSE i najwyższy R². SVR i Elastic Net wykazały niższą dokładność, co może wynikać z nieliniowego charakteru zależności w danych.

Interpretacja modelu za pomocą SHAP

if (!requireNamespace("fastshap", quietly = TRUE)) {
  install.packages("fastshap")
}
library(fastshap)

pred_wrapper <- function(object, newdata) {
  predict(object, newdata = newdata)
}

set.seed(123)
shap_values <- fastshap::explain(
  rf_model,
  X = X_test,
  nsim = 50,
  pred_wrapper = pred_wrapper
)

shap_df <- as.data.frame(shap_values)

shap_summary <- shap_df %>%
  pivot_longer(everything(), names_to = "Feature", values_to = "SHAP") %>%
  mutate(AbsSHAP = abs(SHAP)) %>%
  group_by(Feature) %>%
  summarise(MeanAbsSHAP = mean(AbsSHAP, na.rm = TRUE)) %>%
  arrange(desc(MeanAbsSHAP))

ggplot(shap_summary, aes(x = reorder(Feature, MeanAbsSHAP), y = MeanAbsSHAP)) +
  geom_col(fill = "#d7191c", alpha = 0.9) +
  coord_flip() +
  labs(title = "Ważność cech według SHAP",x = "Cecha",y = "Średnia wartość |SHAP|") + theme_minimal()

shap_summary %>%
  kable(digits = 2, caption = "Tabela 16. Ważność cech według wartości SHAP",
        col.names = c("Cecha", "Średnia SHAP")) %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Tabela 16. Ważność cech według wartości SHAP
Cecha Średnia SHAP
Potential Window (V) 138.22
Upper Limit of Potential Window (V) 55.37
Electrolyte Chemical Formula 54.77
Electrolyte Concentration (M) 54.41
Lower Limit of Potential Window (V) 39.58
Ratio of ID/IG 30.38
SSA 27.22
Current Density (A/g) 17.44
Electrolyte Ionic Conductivity 15.15
Cell_Config 10.54
if (nrow(shap_df) > 0) {
  single_shap <- shap_df[1, ]
  bar_data <- data.frame(
    Feature = names(single_shap), 
    SHAP = as.numeric(single_shap)
  ) %>%
    arrange(desc(abs(SHAP)))
  
  ggplot(bar_data, aes(x = reorder(Feature, SHAP), y = SHAP, fill = SHAP > 0)) +
    geom_col() +
    coord_flip() +
    scale_fill_manual(values = c("TRUE" = "#d7191c", "FALSE" = "#2c7bb6"), labels = c("Zwiększa", "Zmniejsza")) +
    labs(title = "Dekompozycja SHAP dla przykładowej obserwacji", x = "Cecha",  y = "Wpływ SHAP na predykcję",fill = "Kierunek wpływu") +
    theme_minimal()
}

Analiza SHAP potwierdza kluczową rolę parametrów testowych (przedział wartości potencjału, wzór chemiczny elektorlitu, jego stężenie) w predykcji pojemności. Dodatnie wartości SHAP wskazują na wzrost predykcji pojemności, podczas gdy ujemne na jej spadek. Przykładowa dekompozycja pojedynczej obserwacji ilustruje, jak poszczególne cechy przyczyniają się do finalnej predykcji modelu.

Wnioski

Podsumowanie wyników

Zbiór po czyszczeniu zawierał 908 obserwacji i 12 zmiennych bez wartości brakujących, co stanowi 98% pierwotnego zbioru. Wszystkie zmienne numeryczne wykazują istotne odchylenia od rozkładu normalnego, z przewagą rozkładów prawostronnie skośnych i wysoką kurtozą. Zależności między zmiennymi są przeważnie słabe. Najsilniejsze korelacje dotyczą parametrów okna potencjału. Random Forest osiągnął najlepsze wyniki (R² = 0.7), przewyższając SVR i Elastic Net. Analiza SHAP wykazała dominującą rolę zmiennych opisujących warunki testowe (przedział wartości potencjału, wzór chemiczny elektrolitu) nad właściwościami materiałowymi w określeniu pojemności. Model nieliniowy (Random Forest) okazał się znacznie skuteczniejszy niż metody liniowe (SVR, Elastic Net), co potwierdza złożoność zależności między cechami a pojemnością. Wyniki podkreślają kluczowe znaczenie optymalizacji warunków eksperymentalnych (zakres potencjału, typ i stężenie elektrolitu) dla uzyskania wysokiej pojemności elektrod grafenowych.

Ograniczenia analizy:

  • Brak wielu istotnych parametrów materiałowych, które zostały usunięte z powodu zbyt dużej ilości braków w danych
  • Wysoka heterogeniczność danych (różne metody syntezy, warunki testowe)
  • Ograniczona reprezentacja materiałów o wysokiej pojemności (>1000 F/g)